Data exploration

From Lizard guts: Naxos

Kate Culhane, Dept. of Ecology, Evolution, and Marine Biology, UC Santa Barbara
kathrynculhane@ucsb.edu

---

SUMMARY

Exploring the data

---

DEPENDENCIES

##### PACKAGES #####

library(tidyverse) # data manipulation & visualization
library(kableExtra) # pretty tables
library(vegan) # community analyses
library(goeveg) # scree plot for NMDS

##### SOURCE CODE #####

source("R_scripts/source_scales.R") # scales to use in ggplots
source("R_scripts/source_factor_levels.R") # functions to assign factor levels

##### DATA #####

# Abundance data
prey <- read_csv("output/00_prey.csv") %>% # sticky + pitfall samples
  fl_sample()
diet <- read_csv("output/00_diet_inverts.csv") %>% # diet samples
  fl_sample()
inverts <- read_csv("output/00_inverts_all.csv") %>% # ALL invert samples
  fl_sample()

# Site by species matrices
prey_matrix <- read_csv("output/00_prey_matrix.csv") # sticky + pitfall samples
diet_matrix <- read_csv("output/00_diet_matrix.csv") # diet samples
inverts_matrix <- read_csv("output/00_inverts_matrix.csv") # ALL invert samples

Abundance

Number of invert individuals per order

inverts %>% 
  group_by(order, sample_type) %>% 
  summarise(n = sum(count)) %>% 
  pivot_wider(names_from = sample_type, values_from = n) %>% 
  replace(., is.na(.), 0) %>% 
  mutate(trap_total = pitfall + sticky,
         diet = cell_spec(diet,
                          bold = ifelse(diet > 10, "TRUE", "FALSE"),
                          color = ifelse(diet > 10, "red", "black"))) %>% 
  arrange(-trap_total) %>% 
  relocate(diet, .after = last_col()) %>% 
  kbl(escape = FALSE) %>% 
  kable_minimal(full_width = FALSE, position = "left")

order	sticky	pitfall	trap_total	diet
diptera	1808	102	1910	9
hymenoptera	102	1290	1392	129
coleoptera	71	427	498	119
hemiptera	121	148	269	218
acari	9	158	167	14
collembola	1	116	117	0
araneae	15	62	77	76
embioptera	55	10	65	0
raphidioptera	19	0	19	1
thysanura	0	16	16	4
isopoda	0	11	11	9
coleoptera_larvae	0	10	10	0
psocoptera	10	0	10	0
orthoptera	6	3	9	28
chilopoda	0	7	7	0
pseudoscorpiones	0	6	6	2
lepidoptera	3	2	5	10
opilliones	0	4	4	0
lepidoptera_larvae	0	3	3	0
gastropoda	1	1	2	33
homoptera	1	0	1	0
isoptera	0	1	1	0
neuroptera	0	1	1	0
blattodea	0	0	0	5
egg	0	0	0	10
larvae	0	0	0	54
odonata	0	0	0	1

• Number of invert individuals summed across samples
• Orders with >10 individuals in diet samples highlighted in red

Trapping samples

# Abundance by order
prey %>% 
  ggplot(aes(x = sample_type, y = count)) +
  geom_jitter(aes(color = sample_type),
              width = 0.1, size = 1.5) +
  scale_color_manual(values = scol_trap) +
  geom_boxplot(fill = NA, outlier.shape = NA) +
  facet_wrap(~ fct_reorder(order, count, .fun = sum, .desc = TRUE),
             scales = "free_y") +
  theme_classic() +
  theme(legend.position = "NA") +
  labs(x = "Trap type", y = "Abundance (indv/trap)")

• Number of invert individuals trapped per order
• Each point represents a trap (12 trap pairs/site, 9 sites)

Diet samples

# Abundance by order
diet %>% 
  # Data wrangling
  group_by(order) %>% 
  mutate(n = length(count)) %>% 
  # Plot
  ggplot(aes(x = fct_reorder(order, count, .fun = sum, .desc = TRUE), y = count)) +
  geom_jitter(width = 0.1, size = 1.5) +
  geom_boxplot(fill = NA, outlier.shape = NA) +
  geom_text(aes(y = 10, label = paste("n =", n)),
            stat = "unique", angle = 90, size = 3, hjust = 0) +
  scale_y_continuous(expand = c(0,0),
                     limits = c(0,11.5),
                     breaks = c(0,2,4,6,8,10)) +
  theme_classic() +
  theme(legend.position = "NA",
        axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1)) +
  labs(x = "Order", y = "Abundance (indv/lizard)")

• Number of invert individuals consumed per order
• Each point represents a stomach sample from a single lizard (148 lizards)
• Four outliers not plotted for hemiptera: 48, 22, 22, 17

NMDS

Invert trap data

# Wrangle matrix
m_p <- as.matrix(prey_matrix[-c(1:4)])
rownames(m_p) <- prey_matrix$sample
m_p_hell <- decostand(m_p, 'hellinger') # Hellinger transformation
m_p_meta <- select(prey_matrix, sample, site, trap, sample_type) %>% # metadata
  fl_sample()

# Scree plot to check stress per number of dimensions
dimcheckMDS(m_p_hell, distance = "bray", autotransform = FALSE, k = 10)

# Create NMDS ordination
ord_p <- metaMDS(m_p_hell,
                 distance = 'bray', # use Bray-Curtis distances
                 autotransform = FALSE, # already manually transformed the matrix
                 k = 3, # number of dimensions
                 trymax = 1000)

# Stress plot
stressplot(ord_p)

# Plot ordination
par(mar = c(4,4,1,1))
plot(ord_p, display = 'species', type = 'n')
points(ord_p, display = 'sites', pch = 19,
       col = scale_trap[m_p_meta$sample_type])
ordiellipse(ord_p, groups = m_p_meta$sample_type,
            label = TRUE, col = scale_trap, lwd = 2)
text(ord_p, display = 'species')

• Each dot is a trap
• Number of dimensions = ord_p$ndim
• Stress = round(ord_p$stress,3)

Diet + trapping data

# Wrangle matrix
m_i <- as.matrix(inverts_matrix[-c(1:3)])
rownames(m_i) <- inverts_matrix$sample
m_i_hell <- decostand(m_i, 'hellinger') # Hellinger transformation
m_i_meta <- select(inverts_matrix, sample, site, sample_type) # metadata

# Scree plot to check stress per number of dimensions
dimcheckMDS(m_i_hell, distance = "bray", autotransform = FALSE, k = 10)

# Create NMDS ordination
ord_i <- metaMDS(m_i_hell,
                 distance = 'bray', # use Bray-Curtis distances
                 autotransform = FALSE, # already manually transformed the matrix
                 k = 2, # number of dimensions
                 trymax = 1000)

# Stress plot
stressplot(ord_i)

# Plot ordination
par(mar = c(4,4,1,1))
plot(ord_i, display = 'species', type = 'n', xlim = c(-0.15, 0.15))
points(ord_i, display = 'sites', pch = 19,
       col = scol_sample[m_i_meta$sample_type])
ordiellipse(ord_i, groups = m_i_meta$sample_type,
            label = TRUE, col = "black", lwd = 2)
text(ord_i, display = 'species')

• Each dot is a trap
• Number of dimensions = ord_f$ndim
• Stress = round(ord_f$stress,3)

---

SESSION INFO

devtools::session_info()

## - Session info --------------------------------------------------------------------------------------
##  setting  value                       
##  version  R version 4.1.0 (2021-05-18)
##  os       Windows 10 x64              
##  system   i386, mingw32               
##  ui       RStudio                     
##  language (EN)                        
##  collate  English_United States.1252  
##  ctype    English_United States.1252  
##  tz       America/Los_Angeles         
##  date     2021-08-02                  
## 
## - Packages ------------------------------------------------------------------------------------------
##  ! package     * version date       lib source                                 
##  P assertthat    0.2.1   2019-03-21 [?] CRAN (R 4.1.0)                         
##  P audio         0.1-7   2020-03-09 [?] CRAN (R 4.1.0)                         
##  P backports     1.2.1   2020-12-09 [?] CRAN (R 4.1.0)                         
##  P beepr       * 1.3     2018-06-04 [?] CRAN (R 4.1.0)                         
##  P broom         0.7.8   2021-06-24 [?] CRAN (R 4.1.0)                         
##  P cachem        1.0.5   2021-05-15 [?] CRAN (R 4.1.0)                         
##  P callr         3.7.0   2021-04-20 [?] CRAN (R 4.1.0)                         
##  P cellranger    1.1.0   2016-07-27 [?] CRAN (R 4.1.0)                         
##  P cli           3.0.0   2021-06-30 [?] CRAN (R 4.1.0)                         
##  P cluster       2.1.2   2021-04-17 [?] CRAN (R 4.1.0)                         
##  P colorspace    2.0-2   2021-06-24 [?] CRAN (R 4.1.0)                         
##  P crayon        1.4.1   2021-02-08 [?] CRAN (R 4.1.0)                         
##  P DBI           1.1.1   2021-01-15 [?] CRAN (R 4.1.0)                         
##  P dbplyr        2.1.1   2021-04-06 [?] CRAN (R 4.1.0)                         
##  P desc          1.3.0   2021-03-05 [?] CRAN (R 4.1.0)                         
##  P devtools      2.4.2   2021-06-07 [?] CRAN (R 4.1.0)                         
##  P digest        0.6.27  2020-10-24 [?] CRAN (R 4.1.0)                         
##  P dotCall64     1.0-1   2021-02-11 [?] CRAN (R 4.1.0)                         
##  P dplyr       * 1.0.7   2021-06-18 [?] CRAN (R 4.1.0)                         
##  P ellipsis      0.3.2   2021-04-29 [?] CRAN (R 4.1.0)                         
##  P evaluate      0.14    2019-05-28 [?] CRAN (R 4.1.0)                         
##  P fansi         0.5.0   2021-05-25 [?] CRAN (R 4.1.0)                         
##  P farver        2.1.0   2021-02-28 [?] CRAN (R 4.1.0)                         
##  P fastmap       1.1.0   2021-01-25 [?] CRAN (R 4.1.0)                         
##  P fields        12.5    2021-06-25 [?] CRAN (R 4.1.0)                         
##  P forcats     * 0.5.1   2021-01-27 [?] CRAN (R 4.1.0)                         
##  P fs            1.5.0   2020-07-31 [?] CRAN (R 4.1.0)                         
##  P generics      0.1.0   2020-10-31 [?] CRAN (R 4.1.0)                         
##  P ggplot2     * 3.3.5   2021-06-25 [?] CRAN (R 4.1.0)                         
##  P glue          1.4.2   2020-08-27 [?] CRAN (R 4.1.0)                         
##  P goeveg      * 0.5.1   2021-05-10 [?] CRAN (R 4.1.0)                         
##  P gridExtra     2.3     2017-09-09 [?] CRAN (R 4.1.0)                         
##  P gtable        0.3.0   2019-03-25 [?] CRAN (R 4.1.0)                         
##  P haven         2.4.1   2021-04-23 [?] CRAN (R 4.1.0)                         
##  P highr         0.9     2021-04-16 [?] CRAN (R 4.1.0)                         
##  P hms           1.1.0   2021-05-17 [?] CRAN (R 4.1.0)                         
##  P htmltools     0.5.1.1 2021-01-22 [?] CRAN (R 4.1.0)                         
##  P httr          1.4.2   2020-07-20 [?] CRAN (R 4.1.0)                         
##  P jsonlite      1.7.2   2020-12-09 [?] CRAN (R 4.1.0)                         
##  P kableExtra  * 1.3.4   2021-02-20 [?] CRAN (R 4.1.0)                         
##    katereR       0.1.0   2021-07-30 [1] Github (katekathrynkat/katereR@5f2c796)
##  P knitr         1.33    2021-04-24 [?] CRAN (R 4.1.0)                         
##  P labeling      0.4.2   2020-10-20 [?] CRAN (R 4.1.0)                         
##  P lattice     * 0.20-44 2021-05-02 [?] CRAN (R 4.1.0)                         
##  P lifecycle     1.0.0   2021-02-15 [?] CRAN (R 4.1.0)                         
##  P lubridate   * 1.7.10  2021-02-26 [?] CRAN (R 4.1.0)                         
##  P magrittr      2.0.1   2020-11-17 [?] CRAN (R 4.1.0)                         
##  P maps          3.3.0   2018-04-03 [?] CRAN (R 4.1.0)                         
##  P MASS          7.3-54  2021-05-03 [?] CRAN (R 4.1.0)                         
##  P Matrix        1.3-3   2021-05-04 [?] CRAN (R 4.1.0)                         
##  P memoise       2.0.0   2021-01-26 [?] CRAN (R 4.1.0)                         
##  P mgcv          1.8-35  2021-04-18 [?] CRAN (R 4.1.0)                         
##  P modelr        0.1.8   2020-05-19 [?] CRAN (R 4.1.0)                         
##  P munsell       0.5.0   2018-06-12 [?] CRAN (R 4.1.0)                         
##  P nlme          3.1-152 2021-02-04 [?] CRAN (R 4.1.0)                         
##  P permute     * 0.9-5   2019-03-12 [?] CRAN (R 4.1.0)                         
##  P pillar        1.6.1   2021-05-16 [?] CRAN (R 4.1.0)                         
##  P pkgbuild      1.2.0   2020-12-15 [?] CRAN (R 4.1.0)                         
##  P pkgconfig     2.0.3   2019-09-22 [?] CRAN (R 4.1.0)                         
##  P pkgload       1.2.1   2021-04-06 [?] CRAN (R 4.1.0)                         
##  P prettyunits   1.1.1   2020-01-24 [?] CRAN (R 4.1.0)                         
##  P processx      3.5.2   2021-04-30 [?] CRAN (R 4.1.0)                         
##  P ps            1.6.0   2021-02-28 [?] CRAN (R 4.1.0)                         
##  P purrr       * 0.3.4   2020-04-17 [?] CRAN (R 4.1.0)                         
##  P R6            2.5.0   2020-10-28 [?] CRAN (R 4.1.0)                         
##  P Rcpp          1.0.6   2021-01-15 [?] CRAN (R 4.1.0)                         
##  P readr       * 1.4.0   2020-10-05 [?] CRAN (R 4.1.0)                         
##  P readxl        1.3.1   2019-03-13 [?] CRAN (R 4.1.0)                         
##  P remotes       2.4.0   2021-06-02 [?] CRAN (R 4.1.0)                         
##    renv          0.13.2  2021-03-30 [1] CRAN (R 4.1.0)                         
##  P reprex        2.0.0   2021-04-02 [?] CRAN (R 4.1.0)                         
##  P rlang         0.4.11  2021-04-30 [?] CRAN (R 4.1.0)                         
##  P rmarkdown     2.9     2021-06-15 [?] CRAN (R 4.1.0)                         
##  P rprojroot     2.0.2   2020-11-15 [?] CRAN (R 4.1.0)                         
##  P rstudioapi    0.13    2020-11-12 [?] CRAN (R 4.1.0)                         
##  P rvest         1.0.0   2021-03-09 [?] CRAN (R 4.1.0)                         
##  P scales        1.1.1   2020-05-11 [?] CRAN (R 4.1.0)                         
##  P sessioninfo   1.1.1   2018-11-05 [?] CRAN (R 4.1.0)                         
##  P spam          2.7-0   2021-06-25 [?] CRAN (R 4.1.0)                         
##  P stringi       1.6.2   2021-05-17 [?] CRAN (R 4.1.0)                         
##  P stringr     * 1.4.0   2019-02-10 [?] CRAN (R 4.1.0)                         
##  P svglite       2.0.0   2021-02-20 [?] CRAN (R 4.1.0)                         
##  P systemfonts   1.0.2   2021-05-11 [?] CRAN (R 4.1.0)                         
##  P testthat      3.0.4   2021-07-01 [?] CRAN (R 4.1.0)                         
##  P tibble      * 3.1.2   2021-05-16 [?] CRAN (R 4.1.0)                         
##  P tidyr       * 1.1.3   2021-03-03 [?] CRAN (R 4.1.0)                         
##  P tidyselect    1.1.1   2021-04-30 [?] CRAN (R 4.1.0)                         
##  P tidyverse   * 1.3.1   2021-04-15 [?] CRAN (R 4.1.0)                         
##  P usethis       2.0.1   2021-02-10 [?] CRAN (R 4.1.0)                         
##  P utf8          1.2.1   2021-03-12 [?] CRAN (R 4.1.0)                         
##  P vctrs         0.3.8   2021-04-29 [?] CRAN (R 4.1.0)                         
##  P vegan       * 2.5-7   2020-11-28 [?] CRAN (R 4.1.0)                         
##  P viridis       0.6.1   2021-05-11 [?] CRAN (R 4.1.0)                         
##  P viridisLite   0.4.0   2021-04-13 [?] CRAN (R 4.1.0)                         
##  P webshot       0.5.2   2019-11-22 [?] CRAN (R 4.1.0)                         
##  P withr         2.4.2   2021-04-18 [?] CRAN (R 4.1.0)                         
##  P xfun          0.24    2021-06-15 [?] CRAN (R 4.1.0)                         
##  P xml2          1.3.2   2020-04-23 [?] CRAN (R 4.1.0)                         
##  P yaml          2.2.1   2020-02-01 [?] CRAN (R 4.1.0)                         
## 
## [1] C:/Users/kathr/Documents/git-repos/lizard-guts-naxos/renv/library/R-4.1/i386-w64-mingw32
## [2] C:/Users/kathr/AppData/Local/Temp/Rtmp2lbUBn/renv-system-library
## 
##  P -- Loaded and on-disk path mismatch.